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ABSTRACT 

Four turbulence models are described and evaluated for transonic flows using the upwind code CFL3D and 
the central-difference code TLNS3D. In particular, the effects of recent modifications to the half-equation model 
of Johnson-King are explored in detail, and different versions of the model are compared. This model can obtain 
good results for both two-dimensional (2-D) and three-dimensional (3-D) separated flows. The one-equation models 
of Baldwin- Barth and Spalart-AIlmaras perform well for separated airfoil flows, but can predict the shock too far 
forward at the outboard stations of a separated wing. The equilibrium model of Baldwin-Lomax predicts the shock 
location too far aft for both 2-D and 3-D separated flows, as expected. In general, all models perform well for 
attached or mildly separated flows. 


1. Introduction 

As computational fluid dynamics (CFD) grows in 
capability as a tool for the analysis of three-dimensional 
(3-D) aerospace configurations, turbulence modeling 
continues to be one of the primary factors that inhibits 
more widespread usage of Navier-Stokes codes by air- 
craft manufacturers. The search for a model that ac- 
curately predicts both attached and separated 3-D flow 
fields is complicated by the fact that it is difficult to as- 
sess the capabilities of new or refined turbulence mod- 
els because of inherent limitations in the CFD codes 
that use them. Particularly relevant is the issue of trun- 
cation error; the density of the grid used, the type of dif- 
ferencing scheme employed, and, for central-difference 
schemes, the amount of artificial dissipation added for 
numerical stability are all likely to have an effect on 
the solution. 

Turbulence models are often assessed based on 
a limited study - generally only one or two cases 
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are analyzed using a single computer code with little 
or no variation in grid density. As will be shown 
in this paper, such computed results can sometimes 
agree with experimental data almost coincidentally; 
reducing the truncation error through mesh refinement 
or lowering the dissipation levels can sometimes lead 
to dramatically different results, particularly for 3-D 
separated flows. By conducting a broader study that 
covers both two- and three-dimensional configurations 
and using grid refinement studies with more than one 
computer code, the effects of truncation error can be 
determined. Hence, a more accurate assessment of the 
turbulence models is possible. 

This paper investigates some of these issues for 
the studies of transonic flow over the ONERA M6 
wing [1] and the Lockheed Wing C [2]. Two widely 
used computer codes, CFL3D [3] and TLNS3D [4], are 
employed. Each of these codes can employ either the 
Baldwin-Lomax [5] or the Johnson-King [6] turbulence 
model. Additionally, CFL3D can employ the Baldwin- 
Barth [7] or the Spalart-AIlmaras [8] turbulence model. 
To provide a wider basis for assessment, the capabilities 
of the various turbulence models are examined in two 
dimensions using the two-dimensional (2-D) mode of 
CFL3D as well. 
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2. Method 

2.1. The Computer Codes 

CFL3D and TLNS3D both solve the 3-D time- 
dependent thin-layer Navier-Stokes equations with a 
finite-volume formulation. Both can employ grid se- 
quencing, multigrid, and local time-stepping to accel- 
erate convergence to steady state. When converged 
temporally to a steady-state solution, both methods are 
globally second-order accurate. 

CFL3D, described in detail in [3], is an upwind 
code. For all results presented in this paper, upwind- 
biased spatial differencing is used for the inviscid 
terms, and flux limiting is used to obtain smooth solu- 
tions in the vicinity of shock waves. All viscous terms 
are centrally differenced. The equations are solved im- 
plicitly with the use of 3-factor approximate factoriza- 
tion (AF). Either Roe’s flux difference splitting (FDS) 
[9] or van Leer’s flux vector splitting (FVS) [10] can 
be employed to obtain fluxes at the cell faces. 

TLNS3D, described in detail in [4], is a central- 
difference code. Second-order central differences are 
used for all spatial derivatives, and a blend of second- 
difference and fourth-difference artificial dissipation 
terms is used to maintain numerical stability. These 
artificial dissipation terms can be added in either scalar 
or matrix form [11], The solution is advanced explicitly 
in time using either a four-stage or five-stage Runge- 
Kutta time-marching algorithm. 

2.2. The Turbulence Models 

Baldwin-Lomax. The Baldwin-Lomax (B-L) tur- 
bulence model [5] is used widely throughout the CFD 
community; its capabilities and limitations are well- 
known. In short, it is generally considered a good 
model for the prediction of attached flows, but it is de- 
ficient for flows with any significant separated regions. 
In particular, the B-L model tends to predict shocks 
too far downstream for separated transonic flows over 
aerodynamic configurations. 

The B-L model is an algebraic model. The inner 
eddy viscosity is determined via 

Ft#< :P _ L) =p(0.4D ( b-l)J/)^ (1) 

where 

A.b-L) = 1— exp(— »+/26) (2) 

and 

y + = y\fpivTu,/y w (3) 


where y is the distance to the wall or wake cut, £2 is the 
magnitude of the vorticity, and the subscript w refers 
to a wall value. 

The outer eddy viscosity is given by 

Ft,o (B _L) = (0.0168)( 1.6)pT’ ma k e r (B-L) (4) 

where 

F’wake = min [y m F yyi , Z/m(? max (/min )7 F m \ (5) 

and 

r,B— L) = 1/ [l + 5-5(0. 3t//t/ m ) 6 ] (6) 

and 

F = yTlD ( B _L ) (7) 

In wakes, D ( b-L) is taken as 1. The subscript m 
in equations (5) and (6) refers to the value at the 
location where F is at its maximum along grid lines 
that are oriented normal to the body (i.e., in the viscous 
direction for thin-layer Navier-Stokes). The symbol q 
denotes the total velocity. 

The eddy viscosity //, is determined by marching 
away from the wall along these same grid lines and 
is taken as pt )? ; (B _ L) from the wall (where the value is 
zero) to the point above the wall where pt l i {B _ h) first 
exceeds //, 0||l _ l i . Thereafter //, is given the value 
//, 0||l _ l i . Transition is modeled by setting //, to zero 
along all grid lines that are normal to the wall within 
a preselected range. 

Johnson-King. The Johnson-King (J-K) turbu- 
lence model is a so-called one-half equation turbulence 
model that requires the solution of an ordinary differen- 
tial equation (ODE) for the maximum Reynolds shear 
stress over the body surface. Since its introduction 
[6], several modifications and enhancements have been 
made [12], [13], [14], [15]. An attempt will be made 
to describe the different versions of the model, as well 
as the effects of each of the individual modifications. 

In one of the original versions of the J-K model 
[12], (termed J-K1985), the inner and outer eddy vis- 
cosities are given by 


Pt,i { j_ K ) = 0.4p(D ( j_K ))>(!,,„ 

(8) 

and 


Ft, 0(J _ K) = O.OlGSpg^r, j_ K )<7 

(9) 

where 


A J-K) = 1 - exp [-piu'yiiT / (A + y w )] 

(10) 
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and 


r,j-K, = l/[l + 5.5(y/«) 6 ] (11) 

and a ,,, and ut are velocity scales which are both taken 

to be the square root of the maximum Reynolds shear 
2 

stress { — u m ) along each grid line that is normal 
to the body. The displacement thickness <5* is defined 
by 

<5 

S* = (l-q/qt) .(12) 


where 6 represents the boundary-layer edge. A + is 
taken to be 15. 

The term a is the modeling parameter that pro- 
vides the link between the eddy-viscosity distribution 
and the rate equation for the development of the maxi- 
mum Reynolds shear stress. The rate equation is given 
by 


D(-u m ) 0.25, , 

-= (-U 


D 


J" ( 11 in ) ( 11 in ) ~ | 0.25 1) , 

The turbulent diffusion term D m is given by 


(13) 


Dm — 


0.5 (- ,„) " (a 2 - 1) 

0.25(0. 6 - y m ) 

and , n is the dissipation length scale 

0.4t/ m y m /s 0.225 

~~ 0.0 6 y m /6 0.225 


(14) 


(15) 


In equation (13), the subscript eq denotes the equilib- 
rium value of the maximum Reynolds shear stress. In 
equations (14) and (15), y m denotes the distance from 
the body to the location where —u is maximum. 

As discussed in references [6] and [13], the ODE 
(13) can be linearized and solved for , where is 
defined as 


(~U m ) ~ (16) 

The resulting value of is then used to update u 
through an iteration process. The inner and outer 
models are blended to obtain the eddy viscosity //, 
using 


J- K) n 

Pt=Pt,o (J _ K) 1 — exp (17) 

J-K1985 does well in the computation of 2-D 

separated-flow airfoil cases [16], but tends to predict 


shocks too far forward for transonic attached equilib- 
rium flows. Another disadvantage of this version of the 
model is that the determination of the boundary-layer 
edge 6 is computationally difficult and rarely foolproof. 
In fact, this problem is considered by many to be a pri- 
mary disadvantage of the Clauser type of outer eddy- 
viscosity formulation and is one of the reasons why the 
B-L model was developed [5], 

Abid et al. [13] extended the J-K model to 3-D 
flows. Their version (termed J-K1990A) uses equation 
(8) for the inner model, but the velocity scale ut in 
equation (10) is given by 


«t = m | u m , Tw/pku J (18) 

2 

where u m is still the velocity scale (— u m ) . This 

modification to ut is done to provide better predictions 
of skin friction for favorable and zero pressure-gradient 
conditions [14], The J-K1990A model also takes A + 
as 17 and replaces equation (9) with the B-L expression 
(equation (4)), which is modified by the factor cr: 


K = (0.0168)(1.6)/9T 1 ma ker ( B-L)O' (19) 


By using this outer model, Abid et al. avoid the need 
to find 6. The turbulent diffusion term (equation (14)) 
is also modified to contribute to the rate equation only 
in regions of flow recovery (regions where a 1): 


D„ 


0.5 (- ,„) m (a 2 — 1, 0)c 

0.25(0. 6-y m ) 


( 20 ) 


The dissipation length scale is taken to be m = 
min(0.4t/,„, 0.0 6). 

The J-K1990A model was used in reference [13] 
to predict 3-D separated flow over the ONERA M6 
wing successfully on a grid size of 28 65 4 . 

Based on a grid refinement study performed with the 
B-L model, this grid was determined to be fine enough 
to yield grid-converged solutions. However, the fact 
that the model predicted the surface pressures accu- 
rately was largely fortuitous. First, an error in the orig- 
inal coding of the Johnson-King turbulence model in 
TLNS3D moved shocks upstream. Second, computa- 
tions were originally performed with TLNS3D using 
scalar dissipation. A reduction in the dissipation lev- 
els (through the use of matrix dissipation) can have 
a dramatic effect on a separated-flow solution, even 
when the results using scalar dissipation appear to be 
grid converged. Third, recent computations using both 
CFL3D and the corrected version of TLNS3D with ma- 
trix dissipation indicate that the particular ONERA M6 
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case in question may, in fact, be unsteady. Some re- 
sults and further discussion will be given in the results 
section. 

Johnson and Coakley’s modification [14] to the J- 
K model (termed J-K1990J) introduces a blending of 
the inner model between pt,i {3 _ K) (equation (8)) and 
an equilibrium-type model. This blending overcomes 
the tendency for J-K1985 to incorrectly predict shock 
locations for attached equilibrium flows. The blended 
inner eddy viscosity is given by 

Pt,i { = (1 - Pt,i { + /)t,i (J _ K) (21) 

where 

[y T iF 2 + T m 2 / y m T w (22) j 

The velocity scale ir/ is given by equation (18), but the 
velocity scale u m is modified to include a compress- 
ibility correction factor so that 

11 m — Pm/p ( U n , J ~ (23) 

The outer model is the same as in J-K1985 (equation 
(9)); however, a hyperbolic tangent (tanh) blending 
function is employed to merge the inner and outer 
models so that 

Pt=Pt,0 {3 _ K) (pt,i { /Pt,0 {3 _ K) ) II (24) 

This blending is done in place of the exponential (exp) 
type of blending (equation (17)) to give better predic- 
tions of skin friction [14]. The J-K1990J model also 
takes A + as 17 and employs equation (20) for the tur- 
bulent diffusion term. 

The latest published modification to the J-K model 
[15] (termed J-K1992) is primarily a new method for 
determining 6 as a function of F/F, n , or 

6=l.2y m=0 .5, ( (25) 

for the Clauser type of outer model (equation (9)). 
Here, F is the Baldwin-Lomax parameter defined in 
equation (7), and m denotes maximum. With this 
expression, the difficulties that are normally associated 
with computing 6 are avoided. Other than this change, 
the J-K 1992 model is identical to J-K1990J. 

In CFL3D and TLNS3D, both the J-K1990A 
and J-K 1992 versions of the Johnson-King turbulence 
model have been implemented. The Reynolds shear 
stress —u'v' (= t / p) is assumed to be given by 

—u'v' = pt^l/p (26) 


With equation (16), the rate equation (13) is reduced to 
a time-dependent linear equation that is solved using a 
multistage explicit Runge-Kutta time-stepping scheme 
for g, as described in reference [13]. The resulting 
value of g is then combined with the actual value of 
—u'v' m in the flow field to update a via 

a t+At = (27) 

The value of a is limited in practice to lie between 
0.1 and 4. In the J-K1992 model, Pt,i {F<l) in equation 
j (21) is taken as pt,i {B _i ) (equation (1)). In wakes, both 
models revert to the B-L methodology. Transition is 
modeled by setting //, to zero along all of the grid lines 
that are normal to the wall within a preselected range. 

The individual effects of the above-mentioned 
changes and modifications to the J-K model on the 
shock location and extent of separation for 2-D and 3-D 
separated transonic flows (in the authors’ experience) 
are listed in table 1. In general, the effects of differ- 
ent diffusion-term treatments and of the inner-model 
blending are relatively small; however, the effects of 
the other two parameters can be significant, depending 
on the case. As expected, the compressibility correc- 
tion factor \J p n , I j> has little or no effect on transonic 
solutions. The effects listed in table 1 are more pro- 
nounced for more-separated cases; also, 3-D solutions 
tend to be more sensitive to changes in the model than 
2-D solutions. > 

As discussed in reference [15], the Clauser outer 
model has a strong experimental foundation. Baldwin 
and Lomax’s original intent was to emulate this expres- 
sion and avoid the use of 6. However, as shown in ref- 
erence [15], the B-L outer-model expression can over- 
predict the eddy-viscosity levels given by the Clauser 
model by a considerable amount; hence, the B-L outer 
model, used in combination with the J-K inner model 
(as in J-K1990A), predicts shocks further aft than the 
Clauser outer model in combination with the J-K inner 
model (as in J-K1992). 

Changing the blending function between the inner 
and outer models can effect the shock location. Both 
the exponential (equation (17)) and hyperbolic tangent 
(equation (24)) blending functions provide a smooth 
blending; however, the tanh-blended curve lies closer 
to the original two curves near the point where they 
cross than the exp-blended curve. Hence, the maximum 
viscosity values given by the tanh function are slightly 
higher, and the shock position is generally predicted 
further downstream. 

For the treatment of the turbulent diffusion term, 
the maximum (max) function was originally introduced 
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in references [13] and [14] because m has a negligi- 
ble effect in regions where a is less than unity. Also, 
this function generally improved the convergence of 
the solution over computations that employ the abso- 
lute value (abs) function there. Although current com- 
putations do not show negligible effects (the shock po- 
sition is usually altered a small amount), the scheme 
does benefit from improved convergence when m is 
employed only in the recovery regions. 

Finally, the use of inner-model blending (equation 
(21)) causes the shock location to move downstream. 
As an additional significant effect, described in refer- 
ence [14], the predicted skin-friction levels aft of the 
shock are reduced and agree better with experimental 
data and other turbulence models. 


Baldwin-Barth. The Baldwin-Barth (B-B) turbu- 
lence model [7] is a one-equation turbulence model de- 
rived from a simplified form of the - equations. The 
model solves a partial differential equation (PDE) over 
the whole field for the modified turbulent Reynolds 
number : 


2 

t 

(28) 

t 

2 — \/ 2 ^ 


Then 


Pt=P 2 (29) 

The variables , 2 . and 2 are damping terms that 

are described in detail in reference [7]. In this for- 
mulation of the B-B turbulence model, the thin-layer 
assumption has been used for the source term (the last 
term in the PDE). 

This equation is solved implicitly using 3-factor 
AF, with first-order upwind differencing used on the 
advective terms. Local time-stepping is employed to 
accelerate convergence. The transition location is mod- 
eled by phasing out the source term along all of the grid 
lines that are normal to the wall within a preselected 
range. 


Spalart-Allmaras. The Spalart-Allmaras (S-A) 
turbulence model is a one-equation turbulence model 
derived “using empiricism and arguments of dimen- 
sional analysis, Galilean invariance, and selective de- 
pendence on the molecular viscosity” [8], The model 
solves a PDE over the whole field for a working vari- 
able related to the eddy viscosity: 

_ 2 


— t2 

~ t’2 2 t2 ~ 

2 



Then 

P / 

Pt = - ' (31) 

The terms t n , 2. and are described in detail in 
reference [8]. The term denotes the nearest distance 
to any wall. The S-A model is very similar in form 
to the B-B model, although the S-A model includes a 
destruction term that is not present in the other model. 
The PDE is solved using the same implicit method used 
in the B-B model, and transition is modeled in a similar 
fashion as well. 

The B-B and S-A models are relatively new; a 
consensus of opinion has not yet been reached on 
their effectiveness in predicting 2-D and 3-D separated 
aerodynamic flows. Also, from private communication 
with Barth and Spalart, both of the models are being 
revised; hence, any conclusions reached in this paper 
regarding the capabilities of the models are temporary. 

Unlike both the B-L and J-K models, the B-B and 
S-A models lend themselves easily to programming on 
unstructured meshes because no inherent dependency 
on grid structure exists (e.g., there is no need to locate 
maximum values of quantities along grid lines normal 
to the wall). Also, no division into an inner and 
outer model, or into a wall and wake model, is made. 
Because they require the solution of a PDE over the 
whole flow field at each time step, these models are 
more expensive than B-L or J-K; however, the solution 
to the PDE need not be converged fully at each time 
step for iteration toward a steady-state solution. The 
net result is an increase of less than 20% in CPU time 
over the B-L model. 
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3. Results 

3.1. Airfoil Cases 

All airfoil cases are performed using the 2-D mode 
of CFL3D with FDS. The cases chosen are case 9 and 
case 10 from reference [17] for the RAE 2822 air- 
foil. The wind tunnel corrections for case 9 (attached 
flow) are essentially those recommended in reference 

[17] , These corrections are also used by many other re- 
searchers, including the authors of references [14] and 

[18] . The corrections for case 10, taken from reference 
[7], are slightly different than those recommended in 
reference [17], 

The primary grid used is a C mesh with 

177 points on the airfoil, a minimum spacing at the wall 
of 0.0000014 (where represents chord), yielding an 
average + over the body for case 10 of about 0.25. 
The outer boundary extent is approximately 15 , and 
transition is assumed at 3% . 

The conditions used for case 9 are = , 

= , and = . This case has little 

to no separation, for which the J-K1985 model predicts 
the shock location too far upstream [18]. Figures 1(a) 
and (b) show results computed with J-K1992 and three 
other turbulence models compared with experimental 
results. Surface pressure coefficients and upper surface 
skin-friction coefficients (nondimensionalized by edge 
values) for all four models are fairly consistent. B- 
L predicts the shock location the furthest downstream; 
and B-B, the furthest upstream. 

The conditions used for case 10 are = , 

= , and = . Figures 2(a) and 

(b) compare the computed results using the B-L, I- 
K1992, B-B, and S-A models to experimental results. 
The shock location is computed consistently by the J- 
K1992, B-B, and S-A models as slightly downstream 
of the experimental result. The B-L model predicts 
the shock location even further downstream. Skin- 
friction coefficients are predicted consistently by the 
four models upstream of the shock, but significant 
variation occurs downstream. In particular, the B-B 
model indicates that the flow does not reattach aft of 
the shock, unlike the other three models. 

These CFL3D solutions are run using a 3-level W- 
cycle multigrid algorithm. The 2 norm of the residual 
of the equation for density generally converges 3 to 
4 orders of magnitude in 400 cycles, although the lift 
generally reaches its steady-state value within 200 to 
300 cycles. The CPU times required on the Cray-YMP 
computer for 400 cycles to obtain the results shown in 


figure 2 are 214, 230, 249, and 254 CPU seconds for 
the B-L, J-K1992, B-B, and S-A models, respectively. 

The effect of grid refinement on the computed 
results for the J-K1992 model is shown in figure 3. 
The coarser grid consists of every other point on the 
fine grid. The surface pressure predictions do not differ 
significantly, as shown in figure 3(a). In figure 3(b), the 
skin friction in front of the shock decreases as the grid 
is refined, but the change is relatively small between the 
two grids. These results indicate that the grid 

is fine enough to adequately remove most grid-related 
truncation error from the computations. Although not 
shown, the effect of grid refinement on the results with 
the other turbulence models is similar. 

The differences between the results using J-K1992 
and I-K1990A for case 10 are shown in figures 4(a) and 
(b). The latter model predicts the shock slightly further 
downstream, in poorer agreement with the experimental 
results, and predicts skin-friction values that are lower 
in front of the shock and higher aft of the shock. This 
comparison in 2-D is shown here because later results 
in 3-D indicate that the J-K1990A model actually per- 
forms better than J-K1992 for a particular separated 
flow case. The difference between the predicted shock 
locations with the J-K1992 and J-K1990A models is 
not as large as the difference reported in reference [15] 
between the J-K1992 model using the Clauser p. t , ver- 
sus the B-L n t . This disagreement is primarily due 
to the use of the exp blending in J-K1990A (as op- 
posed to tanh blending), which partially counteracts 
the tendency of the B-L outer model to move shocks 
downstream. However, the use of exp blending is also 
the primary cause for the lower skin-friction levels in 
front of the shock using J-K1990A. The difference in 
skin friction behind the shock is attributable primarily 
to the different formulations of the inner model. 

3.2. Wing Cases 

Both CFL3D and TLNS3D are used to compute 
the flow over the ONER A M6 wing [1] and the Lock- 
heed Wing C [2], No corrections to the wind tunnel 
test conditions are employed for the ONERA M6 wing 
cases, as recommended in reference [1], For the Lock- 
heed Wing C case, the corrections of Garriz et al. [19] 
are employed. For most of the ONERA M6 runs, a 
C-O mesh is employed with a minimum 
normal spacing over the wing of 0.000015 and a 
distance from the wing to the outer boundary of at least 
7.95 . For the Lockheed Wing C calculations, a 

C-O mesh is employed with a minimum 
normal spacing over the wing of 0.000001 and a 
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distance from the wing to the outer boundary of at least 
9.11 . For all computations, transition is taken at 

approximately 3% . 

A more extensive set of computations is performed 
for the ONERA M6 wing. These include the use of 
three different test conditions and two different grid 
densities. The first case computed is an attached flow 
case ( = , = , and = 

based on the mean aerodynamic chord). For this 
case, the minimum normal spacing of the grid yields 
an average + of about 4 over the body. CFL3D 
results, using FDS, are shown in figures 5(a) through 
(d) for four different span stations, using four different 
turbulence models. All models give similar results 
that are in good agreement with experimental results. 
Although not shown here, CFF3D with the J-K1990A 
model and TFNS3D with the B-F, J-K1992, and J- 
K1990A models give nearly identical results. 

The first separated-flow case that was attempted 
is the = case, also reported in reference [13]. 
Using the J-K1990A model, the solution was found 
to vary extensively for TFNS3D and CFF3D, depend- 
ing upon the levels of dissipation added and the grid 
density. Figure 6(a) shows TFNS3D solutions at span 
station = qh a grid with 

both scalar and matrix dissipation. (Results at other 
span stations are similar.) The scalar dissipation case 
converges in excellent agreement with experimental re- 
sults. However, lowering the dissipation levels gives 
unsteady, massive separation and moves the shock for- 
ward. Although not shown, increasing the grid density 
to using scalar dissipation causes very 

little change from the scalar case. This 

lack of change by itself would suggest a grid-converged 
solution; however, with the solution using matrix dis- 
sipation in mind, this is clearly not the case because 
solutions using different dissipation levels should tend 
toward the same steady-state answer as the grid is re- 
fined. 

CFF3D shows trends similar to TFNS3D, as 
shown in figure 6(b). Here, the FVS scheme, which 
has more inherent dissipation than FDS, results in a 
steady-state solution with the shock position computed 
in good agreement with experimental results. How- 
ever, the FDS solution is unsteady and massively sep- 
arated, with the shock position considerably upstream. 
These computations indicate that either (1) at these 
conditions this wing is indeed massively separated and 
(probably) unsteady, as shown by the computations 
with low dissipation levels, or (2) a much finer grid 
is needed to find a grid-converged solution. These 
computations also demonstrate how easily one can be 


misled into thinking that a code is giving the “cor- 
rect” answer, when in fact the truncation error has 
not been sufficiently reduced through grid-refinement 
and/or lowering dissipation levels. 

Because of questions regarding the nature of the 
ONERA M6 wing case for = , comparisons are 

made instead at = , for which shock-induced 

separation exists at a lesser extent. The Mach number 
for this case is = . Results using CFF3D 

with FDS and four different turbulence models are 
shown in figures 7(a) and (b) at two representative span 
stations. As expected, the B-F model predicts the shock 
too far downstream in the region of the wing where sep- 
aration exists. The B-B model predicts unsteady, mas- 
sive separation outboard of about = a/id 

a shock location too far forward in comparison with 
experimental data, particularly near the wing tip. Both 
the J-K1992 and S-A models predict similar shock loca- 
tions, upstream of experiment at the outboard stations, 
as well as separation over most of the wing behind the 
shock at these same stations. TFNS3D solutions with 
J-K1992 (not shown) show similar behavior. 

For the Johnson-King turbulence model, better 
3-D results for this case can be obtained using J- 
K1990A. In 2-D, this version predicts shocks somewhat 
downstream of the J-K1992 model, and it does not 
predict skin-friction levels as accurately. (Recall figure 
4.) In 3-D, however, the J-K1990A model shows 
excellent convergence properties and good agreement 
with experimental surface pressures. It is difficult, 
therefore, to ignore this model as a viable alternative. 

Figures 8(a) through (d) show the CFF3D and 
TFNS3D results on two different grids for the = 
case, with the use of the J-K1990A model. As 
would be expected for grid-converged solutions, results 
on the two grids for both computer codes are almost 
the same. Although not shown, TFNS3D results with 
scalar dissipation and CFF3D results with FVS on the 
fine grid also give nearly the same results, which pro- 
vides further evidence that grid convergence has been 
achieved for this case. The shock is generally predicted 
slightly downstream of experimental results. Forward- 
and backward-time particle traces over the wing upper 
surface for the TFNS3D fine-grid solution, shown in 
figure 9, clearly delineate the separated regions; flow 
is from top to bottom. 

Convergence properties for CFF3D and TFNS3D 
for this case on the grid are shown in 

figures 10(a) and (b). Figure 10(a) shows the of 
the 2 norm of the change in density p from one cy- 
cle to the next. (This method is the standard way that 
TFNS3D computes residual; CFF3D output has been 
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altered to obtain the same information.) Both codes 
have a similar convergence history, although TLNS3D 
gives slightly lower residuals at the end of 500 cy- 
cles. Figure 10(b) shows the lift-coefficient conver- 
gence histories. The final for CFL3D is 0.44224, 
while TLNS3D gives 0.44225. CFL3D gets to within 

0.5% of its final lift value in 200 cycles, while TLNS3D 
achieves the same level in 174 cycles. However, 
CFL3D (with FDS) runs faster than TLNS3D (with 
matrix dissipation). On the Cray-YMP computer, 500 
cycles with CFL3D requires 4064 seconds; TLNS3D 
requires 5200 seconds. 

The separated-flow Lockheed Wing C case at the 
nominal conditions = , = , and 

= is run at the corrected conditions 

= , = , and = , as 

recommended in reference [19]. Computations are per- 
formed with CFL3D using several different turbulence 
models, as well as with TLNS3D using the J-K1990A 
model. Results are shown in figures 11(a) and (b) at 
two span stations. At = , all /turbulence mod- 

els give roughly the same results, which are in good 
agreement with experimental results. However, near 
the wing tip, the B-L model predicts the shock too 
far aft, and the B-B model predicts it too far forward. 
All of the other models (including both J-K1992 and 
J-K1990A) yield consistent results that are in reason- 
able agreement with experimental results. TLNS3D 
and CFL3D results with J-K1990A are virtually the 
same. 

4. Conclusions 

Four turbulence models have been described and 
evaluated for transonic 2-D and 3-D flows using the 
computer codes CFL3D and TLNS3D. In particular, 
different versions of the Johnson-King model have been 
described and compared. The following observations 
about all of the turbulence models in general can be 
made: 

1. The Baldwin-Lomax model works well for 
attached flows; however, shocks are predicted too far 
downstream for separated flows. 

2. The Johnson-King model (version J-K1992), 
the Baldwin-Barth model, and the Spalart-Allmaras 
model work well for attached flows and 2-D separated 
flows, but can predict the shock too far upstream for 
some 3-D separated flows. 

3. Version J-K1990A of the Johnson-King model 
generally predicts surface pressures for 3-D attached or 
separated flows very well, but does not work as well 


as version J-K1992 for 2-D flows (particularly in the 
prediction of skin friction). 

Results of this investigation also indicate that ex- 
cessive numerical truncation error can lead to an incor- 
rect evaluation of turbulence models. Specifically, the 
use of scalar dissipation (as opposed to matrix dissipa- 
tion) in a central-difference scheme, or the use of FVS 
(as opposed to FDS) in an upwind scheme, can alter 
the character of 3-D separated-flow solutions. By em- 
ploying both central-difference and upwind computer 
codes to a given problem, in addition to performing 
grid sensitivity studies, this type of uncertainty can be 
minimized. 
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SHOCKS UPSTREAM (MORE SEPARATION) 

SHOCKS DOWNSTREAM (LESS SEPARATION) 

Clauser-type outer model 

Baldwin-Lomax-type outer model 

Exp blending of inner and outer models 

Tanh blending of inner and outer models 

Abs function in diffusion term 

Max function in diffusion term 

J-K inner model 

Blended J-K / equilibrium inner model 


Table 1. Effect of modifications to J-K model. 
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0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 

x/c x/c 

a) surface pressure coefficients b) skin— friction coefficients 

Figure 1. Effect of turbulence model on solution for RAE 2822 airfoil 
case 9, = , = , = , grid, CFL3D. 
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x/c x/c 

a) surface pressure coefficients b) skin— friction coefficients 

Figure 2. Effect of turbulence model on solution for RAE 2822 airfoil 
case 10, = , = , = , grid, CFL3D. 





a) surface pressure coefficients b) skin— friction coefficients 

Figure 3. Effect of grid density on solution for case 10, J-K1992 model, CFL3D. 




Figure 4. 


x/c x/c 

a) surface pressure coefficients b) skin— friction coefficients 

Comparison between J-K1992 and J-K1990A models for case 10, grid, CFL3D. 
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c) 2y/B = 0.90 d) 2y/B = 0.95 

Figure 5. Effect of turbulence model on surface pressure coefficients for ONERA M6 
wing, = , = = , grid, CFL3D. 





x/ c 


x/ c 


a) TLNS3D b) CFL3D 

Figure 6. Effect of dissipation levels on surface pressure coefficients for ONERA M6 wing, = 
= , = , J-K1990A model, grid, = . / 




x/ c 


x/ c 


a) 2y/B = 0.65 b) 2y/B = 0.90 

Figure 7. Effect of turbulence model on surface pressure coefficients for ONERA M6 
wing, = , = = , grid, CFL3D. 
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(F) = 289 x 65 x 49 
(C) = 193 x 49 x 33 



x/ c 


a) 2y/B = 0.65 



x/ c 



0.0 0.2 0.4 0.6 0.8 1.0 

x/ c 


b) 2y/B = 0.80 



0.0 0.2 0.4 0.6 0.8 1.0 

x/ c 


c) 2y/B = 0.90 


d) 2y/B = 0.95 


Figure 8. Effect of grid and computer code on surface pressure coefficients for ONERA 
M6 wing, = , = = , J-K1990A model. 



Figure 9. Particle traces over upper surface of ONERA M6 wing, = , 

= , grid, J-K1990A model, TLNS3D. 




cycles cycles 


a) residual histories b) lift coefficient histories 

Figure 10. Convergence histories for ONERA M6 wing, = , 

= , grid, J-K1990A model. 
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x/c x/c 

a) 2y/B = 0.70 b) 2y/B = 0.90 

Figure 11. Effect of turbulence model on surface pressure coefficients for Lockheed 
Wing C, = , = , = , grid. 
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